* Reset settings and initialize log file
launch, path("share/twomonth")

*-------------------------------------------------------------------------------
* Price and Wasserman (2024), "The Summer Drop in Female Employment"
*
* Description: Show seasonal patterns in the prevalence of non-employment spells
*              encompassing two or more CPS reference weeks.
*-------------------------------------------------------------------------------


* Prepare data and run models
*-------------------------------------------------------------------------------

if "$estimate" != "0" {
	* Load data on adult individuals
	gzuse "$basepath/data/derived/cps_bms_sample.dta.gz", clear

	* Tag whether workers are employed and at work
	gen byte atwork = (emp == 1 & absent == 0)

	* Measure non-working spells lasting exactly 1 CPS month
	gen byte w_n_w = (L.atwork == 1 & atwork == 0 & F.atwork == 1) if !missing(L.emp, emp, F.emp)

	* Measure non-working spells lasting 2+ CPS months
	gen byte w_n_n = (L.atwork == 1 & atwork == 0 & F.atwork == 0) if !missing(L.emp, emp, F.emp)

	* Weeks elapsed for lead month
	gen byte f_weeks = F.weeks

	* Linkability in lead month
	gen byte f_linked_monthly = F.linked_monthly

	* Retain individuals who are observed in three consecutive months
	keep if !missing(L.emp, emp, F.emp)

	* Restrict to those with valid longitudinal links
	keep if linked_monthly == 1 & f_linked_monthly == 1

	* Collapse data
	gcollapse (mean) w_n_w w_n_n (rawsum) wtraked (first) month tmspline* weeks f_weeks [pw = wtraked], by(female tm)
	tsset female tm

	* Loop over sex
	foreach f of numlist 0 1 {
		* Loop over outcome variables
		foreach yvar of varlist w_n_w w_n_n {
			* Run the specification
			quietly ivreg2 `yvar' ib5.month D.tmspline* D.weeks D.f_weeks if female == `f' [aw = wtraked], bw($bandwidth) robust small
			process_estimates, path("twomonth") model("f`f'_`yvar'")
		}
	}
}


* Prepare estimates
*-------------------------------------------------------------------------------

* Load estimates into memory
load_estimates, path("twomonth")

* Prepare coefficient labels for each month
make_coeflabels


* Plot seasonality in two-month transitions by gender
*-------------------------------------------------------------------------------

* Prepare background shading
clear
set obs 14
gen rlow = -1.25
gen rupp = 3.25
gen rval = _n - 0.5

* Plot estimates
foreach yvar in "w_n_w" "w_n_n" {
	if "`yvar'" == "w_n_w" {
		local title "Non-working for just 1 reference week"
	}
	else if "`yvar'" == "w_n_n" {
		local title "Non-working for 2+ reference weeks"
	}

	#delimit ;
	coefplot
		(f1_`yvar', keep(coef*) offset(-.05) recast(connected) label("Women"))
		(f0_`yvar', keep(coef*) offset(+.05) recast(connected) label("Men")),
		coeflabels(`coeflabels')
		title("{bf:`title'}")
		xtitle("")
		xlabel(, alternate)
		ytitle("")
		yscale(range(-1.25 3.25))
		ylabel(-1(1)3, format(%5.1f))
		addplot(
			rarea rupp rlow rval if inrange(rval, 0.5, 5.5), color($ltgs) plotregion(margin(t=0 b=0)) below ||
			scatteri 0 0.5 0 13.5, recast(line) color(black) below)
		vertical
		legend(rows(1) size(*.8) order(4 6) label(4 "Women") label(6 "Men"))
		rescale(100);
	#delimit cr

	tempfile g_`yvar'
	graph save "`g_`yvar''"
}

* Create pdf file
grc1leg "`g_w_n_w'" "`g_w_n_n'", l1title("      Difference relative to May (p.p.)", size(*.8)) imargin(vsmall)
nicepdf "$basepath/output/twomonth.pdf", panels(2) indirect replace

* Close the log file
unlaunch
